#coding: utf-8

xi = []
yi = []
zi = []
ri = []

n = 10

def _Ei(i,x,y,z,r):
	return pow( (x-xi[i]),2 ) + pow( (y-yi[i]),2 ) + pow( (z-zi[2]),2 ) - pow(r,2)
	
def _Fk(k,x,y,z,r):
	res = 0
	for i in range(0,n):
		if i!= k:
			res += pow( ( _Ei(k,x,y,z,r) - _Ei(i,x,y,z,r) ),2 )

	return res
	
def _der_F_xk(k,x,y,z,r):
	res = 0
	for i in range(0,n)
		if i!=k:
			res+=(1-_Ei(i,x,y,z,r))
			
	res*=4*_Ei(k,x,y,z,r)*(x-xi[k])
	
	for j in range(0,n):
		if j!=k:
			res+=4*(_Ei(j,x,y,z,r)-_Ei(k,x,y,z,r))*(x-x[j])
	
	return res
	

def _der_F_yk(k,x,y,z,r):
	res = 0
	for i in range(0,n)
		if i!=k:
			res+=(1-_Ei(i,x,y,z,r))
			
	res*=4*_Ei(k,x,y,z,r)*(y-yi[k])

	for j in range(0,n):
		if j!=k:
			res+=4*(_Ei(j,x,y,z,r)-_Ei(k,x,y,z,r))*(y-y[j])
	
	return res	
	
def _der_F_zk(k,x,y,z,r):
	res = 0
	for i in range(0,n)
		if i!=k:
			res+=(1-_Ei(i,x,y,z,r))
			
	res*=4*_Ei(k,x,y,z,r)*(z-zi[k])
	
	for j in range(0,n):
		if j!=k:
			res+=4*(_Ei(j,x,y,z,r)-_Ei(k,x,y,z,r))*(z-z[j])	
	
	return res	
	
def _der_F_rk(k,x,y,z,r):
	res = 0
	for i in range(0,n)
		if i!=k:
			res+=(1-_Ei(i,x,y,z,r))
			
	res*=(-8*_Ei(k,x,y,z,r)*ri[k])
	
	for j in range(0,n):
		if j!=k:
			res+=-8*(_Ei(j,x,y,z,r)-_Ei(k,x,y,z,r))*r[j]
	
	return res	
	
def Ei(i,v):
	return _Ei(i,v[0],v[1],v[2],v[3])
	
def F(v):
	res = 0
	for k in range(0,n): 
		ind = 4*k
		res +=_Fk(k,v[ind],v[ind+1],v[ind+2],v[ind+3])	
		
	return res
	
	
def der_F(v):
	res = [0]*4*n
	for k in range(0,n):
		ind = 4*k
		res[ind] = _der_F_xk(k,
			
	
def soma_v(v1,v2):
	res = [0]*len(v1)
	for i in range(0,len(v1)):
		res[i] = v1[i]+v2[i]
	return res
		
def menos_v(v1):
	res = [0]*len(v1)
	for i in range(0,len(v1)):
		res[i] = -v1[i]
	return res
	
def interno_v(v1,v2):
	res = 0
	for i in range(0,len(v1)):
		res += v1[i]*v2[i]
		
	return res	
	
def escalar_v(v1,k):
	res = [0]*len(v1)
	for i in range(0,len(v1)):
		res[i] = v1[i]*k
		
	return res	
	
def busca_linear(v,dk):
	sigma = 0.8
	lbda = 1.0
	lim = 10
	
	while( interno_v(der_Fk( soma_v( v,escalar_v(dk,lbda) ) ),dk) < sigma*interno_v(der_Fk(v),dk) ):
		if lim == 0:
			break			
		lbda/=2
		lim-=1
	
	return lbda
		
def main():
	tenho n*4-uplas de variaveis
	logo, o gradiente possui n*4-uplas de componentes tb
	
	